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Abstract 

Light propagation in an infinite uniform turbid medium is treated as a Markov stochastic process 
of photons to provide an intuitive framework for photon migration. The macroscopic physical 
quantities of photon migration are shown to be completely determined by the microscopic statistics 
of the photon propagation direction in direction space and a generalized Poisson distribution of 
the number of scattering events that includes an exponential decay absorption factor. A proper 
diffusion solution is derived with an exact time-dependent central position and half width of photon 
migration in this framework. The diffusion coefficient is found to be absorption-independent. 
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Propagation in a multiple scattering (turbid) medium such as the atmosphere is com- 
monly treated by the theory of radiative transfer (see, for example, Chandrasekhar's classic 
text jI). Recent advances in ultrafast lasers and photon detectors for biomedical imaging 
and diagnostics have revitalized this field0, 0, 0]. The basic equation of radiative transfer is 
the elastic Boltzmann equation, a non-separable integro-differential equation of first order 
for which an exact closed form solution is not known except for the case for isotropic scat- 
terers as far as the authors knowj^. Solutions are often based on truncation of the spherical 
harmonics expansion of the photon distribution function or resort to numerical calculation 
including Monte Carlo simulation Cai et. al.Q recently used truncation of a cumulant 
expansion of the photon distribution function to solve the elastic Boltzmann equation. 

In this Letter, photon migration in an infinite uniform medium is treated as a Markov 
stochastic process. The solution to the elastic Boltzmann equation with a point source 
propagating initially along the positive z-axis from the origin of space and time is interpreted 
as the probability of finding a photon at any specified location, direction and time. 

In the stochastic picture of photon migration in turbid media, photons take a random 
walk in the medium and may get scattered or absorbed according to the scattering coefficient 
fig and the absorption coefficient fia of the medium. A normalized phase function, /(s • 
s'), describes the probability of scattering a photon from direction s to s'. The free path 
(step-size) between consecutive events (either scattering or absorbing) has an exponential 
distribution /iTexp(— /irc^) characterized by the total attenuation fiT = fis + fJ'a- At an 
event, photon scattering takes place with probability Hs/i^'T (the albedo) and absorption 
with probability ^a/ yUT- This picture forms the basis for Monte Carlo simulation of photon 
migration. 

Here we will show that this simple picture of a Markov stochastic process of photons 
can be utilized to compute analytically macroscopic quantities such as the average central 
position and half width of the photon distribution. The microscopic statistics of the photon 
propagation direction in direction space (solely determined by the phase function and the 
incident direction) provides a basis for computing the macroscopic quantities at any specified 
time and position. The bridge between them is a generalized Poisson distribution Pn{t), the 
probability that a photon has endured exactly n scattering events before time t (solely 
determined by the scattering and absorption coefficients of the medium). 

Denote the ith. propagation position, direction and step-size of a photon as x^'\ s^'^ and 



2 




Figure 1: A photon moving along n is scattered to n' with a scattering angle and an azimuthal 
angle (p m. a photon coordinate system xyz whose z-axis coincides with the photon's propagation 
direction prior to scattering. XY Z is the laboratory coordinate system. 

The initial condition is x'^^-' = (0,0,0) for the starting point and s^^^ = Sq = (0,0,1) 
for the incident direction. The laboratory Cartesian components of x^*^ and s*^*^ are x^^ and 
s^^ {a = 1, 2, 3). The photon is incident at time to = 0. For simplicity the speed of light is 
taken as the unit of speed and the mean free path /i^^ as the unit of length. 

The scattering of photons takes a simple form in an orthonormal coordinate system 
( 7 |"^m| ' ^) attached to the moving photon itself where n is the photon's propagation 
direction prior to scattering and m is an arbitrary unit vector not parallel to n[see Fig. 
The distribution of scattering angle 6 E [0, vr] is given by the phase function of the medium 
and the azimuthal angle is uniformly distributed over [0,27r). For one realization of the 
scattering event of angles {9, 0) in the photon coordinate system, the outgoing propagation 
direction n' of the photon will be: 

, m - (n • m)n ,nxm.^., ^ 

n = j j smc^cos^+T r sm sm + ncosy. (11 

|n X m| |n X m| 

The freedom of choice of the unit vector m reflects the arbitrariness of the xy axes of the 
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<Si>=g <S3>=g^ 



1/(1 -g) 



Figure 2: The average photon propagation direction (vector) decreases as g"' where g is the 
anisotropy factor and n is the number of scattering events. 

photon coordinate system. For example, taking m = (0, 0, 1), Eq. (P) gives 
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sin 6 cos + S3 cos 6. 
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Here s^*-* etc are stated in the laboratory coordinate system. 

The ensemble average of the propagation direction over all possible realizations of {6, 0) 
and then over all possible s*^*) in Eq. turns out to be {^8*^'+-'^^^ = {^s*^*)^ {cos 6) because 9 
and (p are independent and (cos0) = (sin0) = 0. By recursion, 



^s(")) = (sW) (cos^)" = %" = z{l - giY 



(3) 



where g = {cos 9) = 1 — (71 is the anisotropy factor [see Fig. Ej. 

Using Eq. (0) and recognizing the symmetry obeyed by the x, y and z components of s^^\ 



the correlations between the propagation directions are given by 
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where (72 = f (sm'^9'^. On the other hand, the correlation between and s^*^ (j > i) can 
be reduced to a correlation of the form of Eq. (jU due to the following observation 



9(sl''s^) (5) 



where p(s^^^\s^^^) means the conditional probabihty that a photon jumps from s*-*) 
to s(^). Eq. © is a resuh of the Markov property of the process p{s''^'\s''^') = 
J ds^^-^^p{s^^^\s^^-^^)p{s^^-^^\s^'^) and the fact J ds^^hfp{s^^^\s^^-'^^) = gs^^~^^ from Eq. ©. 
Combining Eqs. (0} and (0), and using the initial condition of s^^\ we conclude 

{sfs^^) = ^-fg^-^[l + Ul-g,y], J>^ (6) 

where constants /i = /2 = —1 and = 2. Here we see that the auto-correlation of the 
X, y, or z component of the photon propagation direction approaches 1/3, i.e., scattering 
uniformly in all directions, after a sufficient large number of scattering {a = (3 and j = i oo 
), and the cross-correlation between them is zero {a ^ (3). 

The connection between the macroscopic physical quantities about the photon distri- 
bution and the microscopic statistics of the photon propagation direction is made by the 
probability Pn{t) that the photon has taken exactly n scattering events before time t (the 



(n + l)-th event comes at t). We claim p„(t) obeys the generalized Poisson distribution [l4 1 



= M-eM-t) ^ (/i.t)"exp(-/..t) 
n\ n\ 

which is the Poisson distribution of times of scattering with the expected rate occurrence 
of fi~^ multiplied by an exponential decay factor due to absorption. Here we have used 
/i^^ = 1 as the unit of length. This form of Pnit) can be easily verified by recognizing first 
that po{t) = exp(— t) equals the probability that the photon experiences no events within 
time t (and the first event occurs at t); and second that the probability Pn+i{t) is given by 

JQ fiT JO n\ Ht (n + lj! 

(8) 

in which the first event occurred at t' is scattering and followed by n scattering events up 
to but not including time t, that confirms Eq. ((Zj) at n -|- 1 if Eq. (jZj) is valid at n. The total 
probability of finding a photon at time t 

oo 

^p„(t) = exp(-/i„t) (9) 

n=0 

decreases with time due to the annihilation of photons due to absorption. 



The average propagation direction (s(t)) at time t is then, 

(S(t)) = ^^=0 i^r.)Pnit) 

Plug Eqs. © and 13) into Eq. dTUD, we obtain 

(s(t)) = zexp{-fisgit) = zexp{-t/lt). (11) 

Here It = — g) is usually called the transport mean free path which is the random- 

ization distance of the photon propagation direction. 

The first moment of the photon density with respect to position is thus 

(x(t)) = r (s(r)) dr = zk [1 - expi-t/k)] , (12) 
Jo 

revealing that the center of the photon cloud moves along the incident direction for one 
transport mean free path It before it stops [see Fig. 0] . 

The second moment of the photon density is calculated as follows. Assume p(s2, t2|si, ti) 
is the conditional probability that a photon jumps from a propagation direction Si at time 
ti to a propagation direction S2 at time ^2 {h > > 0), the conditional correlation of the 
photon propagation direction subject to the initial condition is given by 
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Denote the number of scattering events encountered by the photon at states (si,ti) and 
(82,^2) as ni and n2 respectively. Here n2 > rii since the photon jumps from (si,ti) to 
(82,^2)- Eq. (fT!?jl can be rewritten as 

{Sp{t2)s^{ti)) = ^— — — . (14) 

Z^n2>ni Pn2~ni 1^2 " HjPm [H) 



After a straightforward calculation by utilizing Eq. (jUj) and ((7j), we obtain 

{Sf3{t2)saiti)) = [1 + exp{-iJsg2ti)] exp [lJ,sgi{ti - t2)] . (15) 
The second moment is then 

[xl{t)) = 2 j'^ dt2 fj dh{s^{t2)s^{h)) . (16) 



The diffusion coefficient is obtained from ((x^) — (x) ) /2t, i.e. 



]_ I t ^ gf2(l - exp(-/i^^it)) _ 1 - exp{-fisg2 t) 
3t 1 fJ^sgi l^l9i{9i - 92) l^l92{9i - 92) 



^ ^^/^ (35^1 - 5'2)[1 - exp(-;U^^it)] ^ 2[1 - exp(-/x,gf2t)] 3[1 - exp(-/i,gfit)]^ 

3t\/i,^i lJ^l9i{9i - 92) l^l92{9i-92) 2fx'^,gl 

after integration. This exact result does not depend on absorption and agrees with our 
previous independently calculated work [Eq. (21) in Ref. 0]. 

The general form of the photon distribution depends on all moments of the distribution. 
However, after a sufficient large number of scattering events have taken place, the photon 
distribution approaches a Gaussian distribution over space according to the central limit 
theorem^]. This asymptotic Gaussian distribution, characterized by its central position 
and half width {2Dt), is then 



x' + y' {z-{z{t))y 



(18) 



where the normalizing factor C{t) = exp{—fiat) owing to Eq. 0. This provides a "proper" 
diffusion solution to radiative transfer, revealing a behavior of light propagation that photons 
migrate with a center that advances in time, and with an ellipsoidal contour that grows and 
changes shape [see Fig. Ej- 

It is also worth mentioning that the absorption coefficient only appears in the generalized 
Poisson distribution Pnit) through an exponential decay factor exp(— /Xat). This exponen- 
tial factor will be canceled in the evaluation of the conditional moments of the photon 
distribution [see Eqs. (fT!?jl and (jlH)]. Hence, the sole role played by absorption is to anni- 
hilate photons and affects neither the shape of the distribution function nor the diffusion 



coefficient. 
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The results, except for the Gaussian photon distribution Eq. ()18p. are exact under the 
sole assumption of a Markov random process of photon migration. The deviation from a 
Poisson distribution of scattering or absorption events can be dealt with by modifying Pn{t)- 
The Markov random process is usually a good description of scattering due to short-range 
forces such as photon migration in turbid media. In situations where interference of light is 
appreciable, the phase of photon, which depends on its full past history, must be considered 



and this is non-Markovian. One well-known example is weak localization of light. 



12] Non- 



Markov processes may also occur in scattering involving long-range forces such as Coulomb 
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Figure 3: The center of a photon cloud approaches It along the incident direction and the diffusion 
coefficient approaches It/ 3 with increase of time. 



interaction between charged particles in which the many-body effect can not be ignored. 
However, the idea presented here may still be helpful. 

In summary, the macroscopic physical property of photon migration in turbid medium 
has its root in the microscopic statistics of photon propagation direction in direction space 
which is solely determined by the phase function. A generalized Poisson distribution function 
determined by the scattering and absorption coefficients of the medium serves as a bridge 
to connect the microscopic statistics to the macroscopic property. This provides us a clear 
and comprehensive physical picture of photon migration in turbid medium. 
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